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Section  1 


INTRODUCTION 


1 . 1 Monte  Carlo  Method 

With  longevity  or  durability  requirements  being  incorporated 
in  the  structural  design  of  aircraft,  analysis  of  structural  fail- 
ure due  to  fatigue  has  become  an  important  design  criterion.  One 
aspect  of  structural  life  analysis  is  the  growth  of  minute  material 
flaws  in  the  stressed  material  which  after  some  time  grow  to  a 
critical  size  due  to  load  cycling  under  normal  flight  conditions. 

This  fatigue  load  history  for  an  aircraft  is  the  sum  of  several 
load  sources.  Ground-air-ground  cycles,  wind  gusts,  maneuver  loads 
and  engine  noise  all  contribute  to  the  fatigue  loading  and  may  cause 
material  flaws  to  grow  to  a dangerously  critical  size. 

Current  methods  of  predicting  when  the  cracks  will  grow  critical, 
or  what  size  they  will  reach  after  a given  number  of  flight  hours, 
leave  much  to  be  desired  since  even  the  model  of  crack  growth  due 
to  simple  constant  sinusoidal  amplitude  fatigue  is  a controversial 
subject.  The  present  scheme  of  crack  growth  prediction  due  to 
expected  flight  loadings  is  to  choose  a crack  growth  model  and  a 
"representative"  load  spectrum  for  one  or  more  flights.  Crack 
growth  due  to  the  pseudorandom  history  is  then  computed  over 
several  flights.  An  example  of  such  a digital  computer  program  is 
the  Air  Force  CRACKS  II  program  [1] . Results  from  such  a scheme 
can  give  accurate  solutions  for  many  deterministic  cases,  but  a 
more  precise  and  truly  random  load  history  is  required  to  properly 
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predict  crack  growth  dispersion.  Such  a program  would  require  the 
random  generation  of  the  load  spectrum  expected  in  real  flight  con- 
ditions containing  ground-air-ground  (GAG)  cycles,  random  or  deter- 
ministic maneuver  loads,  random  gust  loads,  etc.  This  would  involve 
knowledge  of,  for  example,  the  expected  number  or  spectrum  of  maneuver 
loads  as  well  as  a spectrum  model  for  the  turbulent  gust  loads.  Com- 
bined spectra  which  include  all  of  these  sources  may  be  deduced  from 
experimental  data  sampled  from  aircraft  flight  histories. 

Since  each  crack  growth  solution  over  a finite  load  history 
would  be  unique,  many  crack  growths  would  have  to  be  computed  to 
obtain  a distribution  of  crack  size  solutions  which  could  be  used 
in  risk  assessment  evaluation  of  the  structure.  In  addition  to  the 
loads  being  random,  the  sizes  of  the  initial  flaws  or  crack  sizes 
in  the  material  are  also  random  in  nature.  Introduction  of  the 
initial  crack  size  as  an  additional  random  variable  can  also  be 
required  for  more  realistic  growth  calculations.  Other  random 
structural  properties  may  also  be  introduced  into  the  problem,  such 
as  construction  and  inspection  reliability  factors. 

After  the  introduction  of  appropriate  random  variables,  the 
computation  of  crack  growth  due  to  the  simulated  environment  can 
be  repeated  a sufficient  number  of  times  to  represent  the  life  or 
partial  life  of  the  aircraft,  thus  producing  a final  probability 
distribution  of  crack  sizes  in  a particular  structural  member  from 
a Monte  Carlo  simulation.  If  the  computation  of  each  growth  problem 
is  cheap  enough,  the  Monte  Carlo  technique  can  be  used  as  a practical 
approach  to  fatigue  damage  estimation. 
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The  reduced  cost  of  the  analog  approach  to  crack  growth  com- 
putation is  considered  to  be  most  useful  for  trade-off  studies  in 
which  fatigue  damage  risks  may  be  assessed  in  aircraft  fleets  for 
which  usage  (spectrum)  variation  is  considered.  However,  the  limited 
scope  of  the  present  study  required  restriction  to  one  load  spectrum, 
the  primary  objective  having  been  to  demonstrate  the  technical 
feasibility  of  the  analog  method. 

1 . 2 Constant  Sinusoidal  Amplitude  Problem 

It  is  possible  to  obtain  a closed  form  solution  for  the  dis- 
tribution of  crack  sizes  in  a structure  after  constant  amplitude 
fatigue,  given  an  initial  crack  distribution,  by  making  use  of  a 
technique  of  variable  transformation  [2] . Computing  crack  growth 
due  to  the  sinusoidal  fatigue  loading  will  result  in  a shifting  and 
reshaping  of  the  initial  flaw  distribution  as  a function  of  time. 

If  all  cracks  grew  at  the  same  rate,  the  crack  distribution  would 
only  shift  on  the  crack  size  axis  and  keep  its  original  shape,  but 
in  reality,  larger  cracks  grow  faster  than  smaller  ones.  Hence,  a 
reshaping  is  also  seen  in  the  form  of  a wider  dispersion  at  larger 
crack  sizes.  The  variable  transformation  is  accomplished  in  the 
following  way.  Let  aQ  be  a random  initial  crack  size  with  prob- 
ability density  function  fQ(aQ)  and  a = g(aQ),  where  g(aQ)  is  the 

final  size  of  any  initial  crack  of  size  a . 

u o 

Then : 

f(a)  = fQ  tg_1  (a)  ] I or  f (a)  = f q [g-1  (a)  ] | g-2-|_1  (D 


3 


For  example,  if  a Weibull  distribution  for  initial  crack  size  and 


the  Paris  growth  model  [3]  are  assumed  such  that 

<< 


i (\)  = £ e 


where 


a = initial  crack  size  at  N = 0 
o 

a, 3 = given  parameters  for  a Weibull  distribution 
N = number  of  fatigue  cycles 


and 


2r 


2r  r 


cia/dN  = C p (AK)  ~ cf  ^ a ~ c a 


(2) 


(3) 


and  where  and  r are  material  constants,  then  the  Paris  growth 
solution  is  given  by 


s [ a1  r+  (r-l)  CM  ] 


VfL-r) 


a0  - l « — **  (4) 

Now  substituting  Eq.  2 and  Eq.  4 into  Eq.  1 we  have 


^-r  [a1_r+(r-i)CN  j exp 


[a1  r+  (r- 1)  CM] 


-f*) 


(5) 


When  N = 0,  f (a)  reduces  to  f (a  ) . For  N very  large,  f (a)  -*■  0 or 
almost  all  cracks  have  grown  to  critical  size.  See  Fig.  1 for  typ- 
ical solutions.  When  Aa  is  a random  variable,  the  variable  trans- 
formation becomes  quite  difficult  to  use  because  integration  of  the 
growth  model  over  the  random  load  is  tedious.  Thus,  when  consider- 
ing a more  realistic  loading  environment,  the  Monte  Carlo  method 
could  be  employed. 
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The  analog  computer  is  a good  choice  to  facilitate  such  a 
scheme,  provided  that  computation  costs  are  reasonable.  Digital 
computers  will  suffice  for  numerical  integration  of  the  growth 
model  and  random  variable  generation  as  well.  However,  digital 
computer  costs  would  be  quite  high  for  a full-range  parameter  study 
of  initial  crack  sizes,  material  constants  and  a reasonable  range 
of  load  cases.  Since  the  growth  equations  are  primarily  nonlinear 
first  order  differential  equations,  the  high  speed  analog  computer 
offers  a practical  solution  procedure  at  potentially  only  a frac- 
tion of  the  digital  computer  costs.  The  random  load  simulation 
can  be  accomplished  by  introducing  random  noise  generators  with 
appropriate  statistics  for  the  fatigue  environment.  Additional 
generators  can  be  employed  for  other  nonrandom  load  sources.  An 
analog  computer  Monte  Carlo  scheme  has  been  developed  to  calculate 
the  distributions  of  crack  growth  sizes  due  to  random  fatigue  loads. 
The  resulting  computational  costs  are  several  orders  of  magnitude 
less  than  traditional  digital  schemes.  The  development  of  the 
method  is  the  significant  feature  and  not  the  validity  of  any  par- 
ticular crack  growth  model  or  any  specific  random  load  model.  For 
illustrative  purposes,  the  Paris  crack  growth  equation  is  used 
throughout  with  Gaussian  noise  as  the  random  loading  Act. 
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Section  2 


ANALOG  COMPUTER  PROGRAMMING  OF  THE  PARIS  GROWTH  EQUATION 
2 . 1 The  Analog  Computer 

The  development  of  analog  computer  technology  can  be  attrib- 

* 

uted  to  electronics  advances  during  World  War  II  [4]  . Because 

of  the  tremendous  popularity  that  digital  computers  have  acquired 
recently  in  terms  of  data  handling  capacity  and  computation  speed, 
the  analog  has  been  neglected  altogether  in  many  cases  as  an  effi- 
cient computing  device.  Today  there  exist  many  high  speed  repeti- 
tive-operation analog  computers  capable  of  solving  linear  and  non- 
linear higher-order  differential  equations.  Although  digitals  by 
nature  are  capable  of  better  computing  accuracy  than  analog  compu- 
ters, the  analog  is  capable  of  continuous  function  integration  and 
differentiation  which  digital  computers  find  difficult  to  perform. 

In  terms  of  major  computation  costs,  the  analog  computer's 
initial  machine  cost  must  be  compared  with  digital  computation 
rates  and  rentals . The  cost  of  an  analog  computer  may  range  from 
a simple  twenty-dollar  circuit  to  a several-hundred- thousand  dollar 
facility,  depending  on  the  computer's  capabilities.  The  analog 
computer  used  in  this  study  is  an  Electronic  Associates,  Inc., 

TR-48  PACE**. 

* 

Reference  4 contains  an  excellent  analog  computer  method  review. 

* "k 

This  facility  is  of  1960  vintage  and  cost  about  $25,000  at  that 
time . 
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The  fundamental  feature  of  most  analog  computers  is  the  patch 
board  which  contains  an  array  of  components  which  are  connected  by 
wires  in  a specific  sequence  determined  by  the  particular  computer 
program.  The  basic  components  consist  of  integrators,  inverters, 
summers,  multipliers,  logic  circuits,  etc.,  which  determine  the 
overall  capability  of  the  computer.  The  fundamental  circuit  com- 
ponent of  the  analog  computer  is  the  operational  amplifier  (OPAMP) 
which  is  used  in  the  summers,  inverters  and  integrators,  as  well 
as  in  other  devices.  By  connecting  a feedback  resistor  across  the 
OPAMP  and  adding  an  input  resistor,  the  OPAMP  becomes  a summer  of 
all  inputs  into  the  input  resistor  with  gain  determined  by  the  ratio 
of  the  feedback  resistor  to  the  input  resistor.  The  inverter  is 
actually  a summer  with  only  one  input  since  the  OPAMP  output  is  the 
negative  of  its  input.  Integration  is  accomplished  by  replacing  the 

feedback  resistor  with  a feedback  capacitor.  If  x(t)  is  the  integra- 

1 T 

tor  input,  the  output  will  be  / x(t)  dt  + Vq  where 

Vq  = initial  condition  of  the  integrator  = x(o) 

R = input  resistance 
C = feedback  capacitance 

If  RC  is  variable,  the  integration  rate  or  time  scale  can  be  adjusted 
such  that  a problem  taking  several  minutes  to  solve  in  real  time 
requires  only  a few  milliseconds  with  proper  RC  values  (time  scaling) . 
Table  1 summarizes  some  of  the  most  common  component  devices  used  in 
analog  programming.  Nonlinear  devices  such  as  multipliers,  squares, 
dividers  and  function  generators  are  also  available  in  most  analog 
component  boards.  The  most  common  nonlinear  circuits  are  listed  in 
Table  2. 
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In  analog  programming,  variable  scaling  is  required  for  accu- 
rate problem  solutions.  The  independent  variables  of  the  equations 
are  represented  on  the  analog  computer  by  time.  The  dependent  vari- 
ables and  their  derivatives  with  respect  to  time  are  represented  by 

voltages.  These  signal  voltages  must  not  exceed  the  maximum  allowed 

* 

machine  voltage  , nor  should  they  change  rapidly  enough  to  exceed 
the  frequency  limitations  of  the  computer  or  recording  equipment. 

In  addition,  the  voltages  should  not  be  as  small  as  the  order  of 
circuit  errors,  and  the  problem  solution  time  should  not  be  so  long 
as  to  waste  computer  efficiency.  By  scaling  the  problem  properly, 
rapid , accurate  solutions  can  be  obtained  by  appropriate  voltage 
variables  and  solution  speeds. 

2 . 2 Paris  Equation  for  Constant  Load 

The  first  step  in  setting  up  the  Monte  Carlo  method  is  to 
choose  a crack  growth  model.  The  Paris  equation  has  been  chosen 
for  the  sake  of  simplicity.  The  Paris  growth  equation  is  given  by: 

nI  p Uf  /Z. 

Ja/JH  = C?(AK)  r * Air  a (6) 

Typical  values  for  N and  C for  aluminum  are  N - 4 and  C = 

p p P P 

-10  -9 

3.49  x 10  to  3.91  x 10  , where  a is  in  inches  and  A a is  in 

* * 

ksi  . N has  been  chosen  as  a multiple  of  2 so  that  circuit 
squaring  devices  can  be  conveniently  used  in  the  problem  patching. 


* 

In  the  case  of  the  TR-48,  +10  volts. 

k k 

Parameters  obtained  by  fitting  experimental  da/dN  data  [5] . 
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For  a constant  A a,  the  problem  is  quite  trivial  resulting  in 
the  solution: 

11  4 

- = =■ C Aa  N (7) 

a a p 

o r 

where  a is  the  initial  crack  size  at  N = 0.  The  nature  of  this 
o 

family  of  solutions  is  that  of  an  increasing  exponential  curve 
which  presents  three  major  problems  on  the  analog  computer: 

1.  Without  careful  selection  of  the  solution  time  or 
length,  several  orders  of  magnitude  of  crack  growth 
can  result  which  would  obviously  exceed  the  maximum 
voltage  limits. 

2.  Minute  errors  during  the  first  portion  of  the  solu- 
tion can  grow  significantly  by  the  end  of  the  problem. 

4 

3.  Since  the  growth  rate  is  proportional  to  Aa  , very 
rapid  crack  growths  may  result. 

With  proper  time  and  variable  scaling,  all  of  these  difficulties 
can  be  overcome . 

The  symbolic  analog  programming  diagram  for  the  Paris  growth 

equation  for  N = 4 is  given  in  Fig.  2.  A typical  analog  growth 
P 

solution  recorded  on  chart  paper  is  given  in  Fig.  3.  The  values 
of  Aa  and  a are  recorded  on  the  same  two-channel  recorder  for  com- 
parison. A quick  calculation  shows  the  exact  analytical  solution 

for  a crack  to  grow  from  0.1  inch  to  1.0  inch  for  C = 3.488  x 10 

P 

Np  = 4 and  Aa  = 50  ksi  is  4128  fatigue  cycles.  The  analog  solution 
gives  4130,  verifying  the  accuracy  of  the  analog  program.  It  is 
quite  easy  to  see  that  Aa  can  now  be  a random  variable  of  time  rather 
than  just  a constant. 
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2 . 3 Random  Load  Simulation 

One  of  the  advantages  of  using  the  analog  rather  than  the 
digital  computer  for  growth  calculations  is  the  ease  with  which 
random  signals  can  be  generated.  Band-limited  Gaussian  white  noise 
can  be  easily  generated  by  amplifying  the  potential  across  a Zener 
diode.  Ideally,  white  noise  has  a power  spectrum  which  is  constant 
for  all  frequencies.  The  diode  power  spectrum  in  reality  is  band- 
width-limited and  thus  the  random  noise  can  be  thought  of  as  band- 
width-limited white  noise.  Fig.  4 illustrates  a home-built  circuit 
used  to  obtain  Gaussian  noise  for  load  simulation. 

2 . 4 Track-Hold  Noise 

Because  the  Zener  diode  noise  contains  very  high  frequencies, 
the  OPAMPs  in  the  analog  computer  register  overloading  during  inte- 
gration. A technique  used  to  remedy  this  problem  is  to  first  pro- 
cess the  noise  through  a track-hold  circuit.  The  track-hold  cir- 
cuit is  designed  such  that  the  noise  is  tracked  when  a control  volt- 
age is  positive, and  the  last  sampled  noise  voltage  is  held  constant 
when  the  control  voltage  is  negative.  By  allowing  the  control  volt- 
age to  periodically  jump  from  positive  to  negative,  a track-hold 
circuit  is  obtained  which  can  hold  the  noise  for  a short  time, 
sample  the  signal,  hold  the  new  value  and  so  on.  The  noise  sampling 
frequency  is  determined  by  the  control  voltage  frequency.  The  result- 
ing track-hold  signal  is  a series  of  randomly  occurring,  short  straight 
lines  of  uniform  length  (Fig.  5) . The  track-hold  signal  was  tested 
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to  insure  that  the  Gaussian  distribution  of  the  noise  generator 
was  not  altered  by  the  track-hold  circuit. 

For  the  present,  the  Gaussian  signal  x(t)  from  the  Zener  diode 
will  be  called  Aa  (t)  for  illustrative  purposes.  The  structural  stress 
due  to  gust  loading  a(t)  can  be  fitted  to  a Gaussian  distribution, 
but  that  is  not  the  same  as  Aa(t),  which  is  the  peak-to-peak  ampli- 
tude difference  for  a fatigue  cycle.  The  essence  of  the  Monte  Carlo 
method  is  that  its  development  is  not  dependent  on  the  assumed 
Aa(t)  distribution,  which  is  taken  to  be  Gaussian  for  convenience. 

It  is  convenient  to  let  the  computer  determine  the  mean  and 
standard  deviation  of  the  random  load  as  a check  to  insure  that  the 
proper  load  distribution  has  been  modeled.  Since  the  loading  signal 
is  known  to  be  Gaussian,  all  that  is  necessary  to  describe  the  dis- 
tribution is  the  mean  and  first  moment.  For  continuous  functions, 
the  mean  and  variance  are  given  approximately  by: 


The  symbolic  program  for  these  circuits  is  shown  in  Fig.  6.  Com- 
puting the  mean  and  standard  deviation  over  a finite  time  interval 
is  only  valid  when  the  Aa  (t)  is  an  ergodic  stationary  random  pro- 
cess. Ergodic  stationary  processes  are  those  for  which  each 
ensemble  record  is  statistically  equivalent  to  all  other  records 
and  for  which  ensemble  averages  can  be  replaced  by  time  averages 
of  a single  record.  The  limit  of  the  time  average  over  a finite 
interval  from  time  0 to  T,  as  T tends  to  infinity,  will  be  equivalent 
to  the  ensemble  average  for  ergodic  processes. 
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The  track-hold  circuit  serves  two  purposes.  First,  periodic 
sampling  and  holding  of  the  continuous  Gaussian  noise  is  needed  to 
reduce  the  high  overload  frequencies  to  a level  which  the  analog 
circuitry  can  tolerate.  Second,  the  resulting  blocks  of  A o make 
digital  program  verification  convenient  with  regard  to  the  load 
description. 

Inaccuracies  of  the  crack  size  distribution  will  result 
from  insufficient  sampling  of  the  loading.  For  more  realistic 
load  modeling,  the  track-hold  circuitry  would  not  be  used  and 
either  more  expensive  analog  circuitry  with  wider  frequency 
characteristics  dr  frequency  filtering  of  the  signal  would  be 
necessary. 

2 . 5 Paris  Solution  for  Random  Loading 

Solving  the  Paris  growth  equation  with  random  loading  becomes 
a simple  matter  once  the  appropriate  material  constants  and  load 
distribution  statistics  are  determined.  A symbolic  diagram  of  the 
analog  program  appears  in  Fig.  7.  The  solution  of  the  growth  equa- 
tion actually  takes  very  few  components  compared  to  the  load  genera- 
tion and  statistical  circuits.  Fig.  8 illustrates  a growth  solution 
in  real  time  for  an  initial  crack  of  0.1  inch  subjected  to  random 
loading  with  a mean  of  37  ksi  and  standard  deviation  of  3.9  ksi  for 
15,000  cycles.  The  time  average  and  standard  deviation  of  the  loading 
are  computed  over  a 10  second  interval.  With  the  load  history 
recorded  and  inserted  into  the  digital  CRACKS  II  program,  a solution 
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which  can  be  considered  exact  verifies  the  accuracy  of  the  analog 
solution.  This  check  is  done  to  insure  that  the  rapid  changes  in 
the  loading  during  the  track-hold  sample  mode  are  not  disturbing 
the  analog  circuitry.  It  is  now  apparent  that  the  time  scaling 
can  be  altered  such  that  much  faster  computing  is  available,  pro- 
vided that  solution  accuracy  is  retained. 

2 . 6 Repetitive  Operation  (REPOP) 

High  speed  analog  computers  have  the  capability  of  repeating 
a problem  solution  by  use  of  the  repetitive  operation  (REPOP)  mode. 
This  condition  provides  a means  of  switching  the  integrators  between 
the  reset  (initial  condition)  and  the  operate  (integration)  modes. 
The  TR-48  10-mfd  integrator  feedback  capacitors  are  replaced  with 
0.02-mfd  capacitors  to  change  the  problem  time  scale  by  a factor 
of  500.  Normally  the  REPOP  mode  is  used  primarily  for  sweeping 
oscilloscope  output.  However,  for  solutions  in  the  400  to  500 
millisecond  range,  the  REPOP  mode  can  be  used  with  chart  recording 
equipment.  The  REPOP  mode  permits  the  calculation  of  several 
hundred  problems  in  a very  short  time. 

With  such  a decrease  in  solution  time,  there  exist  cases  where 

crack  growth  or  response  is  so  rapid  that  errors  are  induced  in  the 

OPAMPs  as  well  as  in  the  chart  recorder.  A sample  problem  was  run 

as  before,  but  in  the  REPOP  mode  for  an  accuracy  check.  Fig.  9 

illustrates  the  output  from  the  problem  solution.  By  recording  the 

load  history  and  repeating  the  problem  on  CRA.CKS  II,  two  solutions 

were  obtained,  one  for  a = 0.110  inch  and  one  for  a = 0.112  inch. 

o o 
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Since  it  is  difficult  to  read  the  analog  aQ  to  a third  significant 
figure,  the  two  digital  solutions  were  computed  and  found  to  bound 
the  analog  solution  within  the  output  tolerances.  Note  how  the 
solution  accuracy  deteriorates  with  a relatively  small  error  in  the 
initial  crack  size.  The  exponential  nature  of  the  growth  solutions 
demands  careful  attention  when  using  analog  circuitry,  since  over- 
loading may  occur  without  much  disturbance.  The  above  REPOP  solu- 
tion simulated  crack  growth  over  20,000  fatigue  cycles  lasting  0.40 
seconds,  or  with  a time  scaling  of  1 second  equal  to  50,000  fatigue 
cycles.  For  a series  of  200  problems  assuming  20,000  cycles  each 
and  0.50  seconds  for  each  problem  including  reset  time,  the  entire 
run  would  last  100  seconds. 

2 . 7 Errors  and  Special  Circuitry 

Due  to  the  nature  of  the  family  of  crack  growth  solutions, 
several  sources  of  error  are  present.  Perhaps  the  most  signifi- 
cant is  the  potential  range  and  frequency  of  the  exponential  solu- 
tions occurring  in  less  than  one-half  second  of  computation  time. 

4 

Because  the  growth  rates  are  functions  of  the  loading  A a (t) , there 
exist  cases  where  the  loads  will  become  large  enough  during  the 
problem  solution  time  for  crack  sizes  to  grow  rapidly  by  several 
orders  of  magnitude.  This  presents  a difficult  problem  for  analog 
computers  since  their  variable  range  is  only  one  order  of  magnitude, 
or,  at  most,  two.  However,  it  is  possible  to  reduce  this  problem 
by  rescaling  variables  internally  in  the  program.  For  example,  in 
problems  where  initial  flaw  sizes  are  0 [0.005  inch]  and  the  solution 
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time  and  loads  are  such  that  cracks  could  grow  to  0 [ 0 . 0 7 inch],  the 
growth  equation  solutions  can  be  divided  into  two  problems: 

1.  A variable  scaling  of  1.0  machine  volt  equal  to 
0.001  inch. 

2.  A rescaling  of  1.0  machine  volt  equal  to  0.01 
inch  when  the  crack  grows  to  0.01  inch  or  larger. 

If  the  entire  two-magnitude  range  is  attempted  in  the  scaling  of 

1 volt  equal  to  0.01  inch,  the  initial  flaw  value  must  be  less  than 

1 volt.  Beginning  the  solution  at  such  a small  voltage,  especially 

for  an  exponential  solution,  taxes  the  tolerance  limits  of  the  analog 

accuracy  over  the  10-volt  range.  Observing  the  solutions  on  an 

oscilloscope,  it  is  quite  easy  to  see  the  accumulation  of  errors 

at  an  initial  crack  voltage  level  of,  say,  0.5  volts  by  observing 

apparent  growth  for  no  loading.  More  modern  and  sophisticated 

analog  computers  can  have  a variable  range  of  +100  volts,  which 

greatly  reduces  the  problems  of  two-magnitude  solution  ranges. 

However,  for  this  study,  the  TR-48  has  a +10  volt  range  and  special 

circuitry  or  programming  is  required. 

The  design  of  circuitry  to  consider  both  small  and  large  cracks 
is  based  on  the  availability  of  logic  flip-flops  or  electronic 
switches  controlled  by  signal  comparators.  The  flip-flops  are 
solenoid-controlled  switches  which  have  a propagation  time  of 
several  milliseconds  and  should  be  avoided  if  possible,  since  they 
are  not  entirely  reliable.  Electronic  switching  devices  are  much 
more  expensive  but  are  worth  the  extra  cost  since  they  are  quite 
fast  (propagation  time  approximately  5 microseconds) . A logic  cir- 
cuit can  be  devised  with  electronic  switches  such  that  when  the 
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crack  size  is  below,  say,  0.01  inch, the  first  growth  circuit  for 
small  cracks  is  integrating  and  the  second  circuit  for  large  cracks 
is  in  stand-by  in  the  reset  mode  until  the  crack  grows  to  0.01  inch. 
Then,  with  appropriate  switching,  the  first  circuit  is  taken  out  to 
prevent  overloading  and  is  replaced  by  the  rescaled  second  circuit 
with  an  initial  condition  of  0.01  inch  and  allowed  to  continue  inte- 
grating the  solution  for  the  remaining  computing  time.  A third  stage 
can  be  added  if  required,  or  an  overload  prevention  switch  can  be 
used  to  prevent  growth  beyond  10  volts  (0.1  inch)  in  the  new  rescaled 
circuit. 

In  the  REPOP  Monte  Carlo  program,  caution  must  be  taken  to 
prevent  OPAMP  overloading  or  saturation  with  voltages  greater  than 
10  volts,  since  a few  seconds  are  required  to  unsaturate  the  ampli- 
fiers. When  using  loads  with  a large  variance,  many  solutions  might 
overload  the  circuitry  requiring  shut-down  and  relaxation  time  before 
continuing. 

In  addition  to  variable  errors,  inaccuracies  in  REPOP  comput- 
ing time  are  also  significant.  Figure  10  shows  the  REPOP  solution 
of  a typical  problem  in  which  A a is  constant.  The  chart  speed  has 
been  reduced  to  show  several  problems.  Since  Aa  is  constant,  the 
computer  is  solving  the  same  problem  repeatedly,  but  the  final 
solution  values  are  not  consistent.  Due  to  slight  variations  in 
the  REPOP  resetting  and  operating  switches,  the  exponential  solu- 
tions exhibit  variance  due  to  microsecond  propagation-time  errors 
in  the  switches.  Proper  care  must  be  taken  to  insure  constant  solu- 
tion time.  When  chart  recording  equipment  is  used,  each  solution 
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must  be  recorded  by  hand.  In  the  case  of  digital  sampling  equip- 
ment, care  is  needed  to  make  sure  that  identical  solutions  are 
recorded  for  constant  loading. 

2 . 8 Sample  Monte  Carlo  Run 

With  the  program  used  in  Fig.  11  patched  into  the  analog  com- 
puter and  a constant  A a test  solution  made  to  insure  correct  material 
constant  settings,  a sample  Monte  Carlo  run  was  made  for  250  problems. 
The  initial  crack  size  remained  constant  at  0.1  inch  with  a maximum 
allowed  solution  of  1.0  inch.  Using  typical  Paris  constants  for 
aluminum  and  Gaussian  loading  with  a 26.7  ksi  mean  and  7.5  ksi  .stan- 
dard deviation,  crack  growth  distributions  were  recorded  at  10,000, 
15,000,  and  20,000  fatigue  cycles  as  given  in  Figs.  12  through  14. 
Because  the  range  of  solutions  was  well  within  one  magnitude,  no 
special  switching  and  rescaling  circuitry  was  required.  As  expected, 
the  distributions  shift  and  spread  out  with  increased  fatigue. 
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Section  3 


EXTENDED  MONTE  CARLO  SIMULATION 

An  extended  Monte  Carlo  simulation  of  the  Paris  equation  was 
conducted  to  generate  an  example  data  base  for  final  crack  size, 
in  order  to  illustrate  the  use  of  the  Monte  Carlo  method  in  risk 
analysis.  Typical  properties  for  0.09-inch  thick  7075-T6  aluminum 
sheet  [5]  were  used  to  define  the  Paris  equation  constants  and  a 
critical  crack  size.  The  following  Gaussian  distribution  for  stress 
range  A a was  assumed: 


i r 4 / ACT-  1 4.2S~  \ 2-  ~) 

- Tz e*F }'  * V rzr  > j (9) 

where  Aa  is  in  units  of  ksi.  The  simulation  was  repeated  for 
eleven  initial  crack  sizes  ranging  from  0.005  to  0.07  inch.  The 
objective  of  the  simulation  was  defined  to  be  production  of  a 
final  crack  size  histogram,  corresponding  to  each  assumed  initial 
crack  size,  after  growth  computation  over  20,000  load  cycles,  which 
simulated  approximately  4,000  flight  hours.  In  addition,  it  was 
hoped  that  this  amount  of  growth  would  result  in  significant  prob- 
abilities for  achieving  final  crack  sizes  greater  than  the  critical 
size.  The  simulation  achieved  the  latter  objective  partially. 


3 . 1 The  Computation  Problem  and  Its  Solution 


The  Paris  equation  for  crack  growth  may  be  integrated  to  the 


form: 


l-  n p/z>  1-  Np/Z- 

a;  - 


(10) 
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where  is  any  given  initial  crack  size,  and  where  a^  represents 
the  corresponding  final  crack  size  after  N fatigue  cycles.  If  a 
risk  analysis  of  initial  flaw  size  variability  is  intended,  a^  is 
the  independent  random  variable  in  Eq.  10,  and  the  distribution  for 
can  be  derived  in  the  following  simple  manner.  Since  the  load 
history  A a is  deterministic  for  the  purpose  of  the  analysis,  Eq.  10 
represents  a one-to-one  relationship  between  a^  and  a^  values. 
Furthermore,  the  right-hand  side  of  Eq.  10  is  positive  so  that  for 
any  two  a^  values: 

af  > af  if  and  only  if  a.  > a.  (11) 

r2  rl  x2  11 

In  other  words,  the  solutions  do  not  cross  over  as  N is  varied,  and 
therefore : 


Ff(af|ai,  N)  = (12) 

where 

Ff(a^|a^,  N)  = Pr{Final  Crack  Size  £ a^,  given  a^,  N}  (13) 

F^(a^)  = Initial  Crack  Size  <_  a^}  (14) 

Hence,  the  final  crack  size  distribution  may  be  determined  with 
reasonable  accuracy  by  repeating  the  growth  computation  (Eq.  10) 
for  10  to  20  a^  values. 

The  above  direct  procedure  is  made  possible  by  the  explicit 
relation  between  the  dependent  and  independent  random  variables 
a^  and  a^.  However,  the  implicit  relationship  between  a^  and  ho 
prevents  the  inference  of  a simple  analytical  correspondence  similar 
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to  Eqs . 11  and  12.  Also,  an  appeal  to  the  central  limit  theorem  of 
probability  theory  offers  only  a crude  approximation  by  considering 
a^  as  the  sum  of  a large  number  of  uncorrelated  random  variables. 

It  is  therefore  necessary  to  conduct  a Monte  Carlo  experiment  to 
determine  F^  where  service  loading  is  the  independent  random  vari- 
able. If  no  a priori  assumptions  are  made  about  the  form  of  the 
final  crack  size  distribution,  200  to  300  replicate  computations 
are  required  to  obtain  a good  experimental  histogram  [2],  and  2,000 
to  3,000  computations  are  required  to  study  10  selected  a^  values. 

This  heavy  digital  computing  burden  is  increased  further  if  one 
attempts  to  have  the  load  history  simulate  reality  more  closely 
by  increased  fineness  in  layering,  i.e.,  by  assigning  fewer  load 
cycles  per  block  and  increasing  the  randomization  of  block  sequence. 

Analog  computation  was  therefore  examined  and  developed  as 
an  alternate  scheme  for  crack  growth  calculation  in  the  present  pro- 
gram. A brief  review  of  the  analog  scaling  requirements  as  they 
pertain  to  the  relations  between  computing  time,  number  of  events, 
and  equivalent  flight  hours  is  given  here  to  provide  a complete 
picture  of  the  data  base  generation  technique.  The  following  factors 
constitute  of  set  of  constraints  on  the  Monte  Carlo  experiment,  within 
which  a realistic  operating  point  must  be  found. 

First,  efficient  use  of  the  analog  computer  for  several  hundred 
or  more  runs  requires  employment  of  the  repetitive  operation  mode, 
in  which  each  run  consumes  on  the  order  of  0.5  second  in  real  time. 

Of  this  total,  roughly  0.4  second  is  spent  on  computation  and  0.1 
second  is  required  as  reset  time  for  the  next  run.  The  real  time 
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consumed  for  computation  should  correspond  to  the  total  number  of 
load  (Ag)  events  expected  in  one  service  life. 

Second,  the  load  history  is  simulated  by  a random  noise  genera- 
tor and  track-hold  circuit  which  produce  a nearly  square-wave  signal 
for  Ag.  Each  level  portion  of  the  signal  is  of  constant  duration 
(controlled  by  the  track-hold  frequency)  and  may  be  identified  as 
a load  block  representing  a number  of  fatigue  cycles  at  constant 
stress  range.  The  achievement  of  a realistic  simulation  requires 
some  care  in  the  choice  of  track-hold  frequency.  At  one  extreme, 
a very  wide  dispersion  in  af  values  results  if  the  track-hold  fre- 
quency is  very  low.  This  extreme  tends  toward  an  experiment  in 
which  constant-amplitude  fatigue  is  simulated,  with  the  choice  of 
range  being  a random  variable  at  the  beginning  of  each  run. 

This  is  obviously  completely  unrepresentative  of  fleet  service 
experience.  At  the  other  extreme,  the  random  noise  signal  may  be 
input  directly,  without  benefit  of  the  track-hold  filter.  This  is 
desirable  from  a simulation  viewpoint,  but  the  high-frequency  com- 
ponents of  the  signal  tend  to  overload  the  analog  integrating  cir- 
cuits, thus  making  the  results  of  questionable  validity.  Hence, 
a minimum  track-hold  frequency  must  be  chosen,  such  that  at  least 
a few  hundred  samples  of  the  random  noise  are  obtained  during  each 
run  to  insure  adequate  sampling  of  the  stress  range  distribution. 

Combination  of  the  above  requirements  defines  the  relationships 
between  computing  time,  service  life  and  load  blocking.  In  the 
present  case,  the  time  scaling: 

N - 50,000  t (15) 
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was  chosen.  According  to  Eq.  15,  a 0.4-second  computing  time  is 
equivalent  to  20,000  events  or  approximately  3,560  flight  hours. 
This  time  scaling  was  chosen  to  achieve  a reasonable  simulation  of 
one  airplane  service  life. 

Preliminary  experiments  with  track-hold  frequency  indicated 
that  the  dispersion  of  final  crack  size  values  changed  very  little 
for  frequencies  between  200  and  1,200  Hz.  The  track-hold  frequency 
was  set  at  500  Hz  for  the  Monte  Carlo  experiment;  i.e., 

500  Hz  x 0.4  sec  = 200  samples  per  run  (16) 


50,000  events/sec 
500  Hz 


= 100  events  per  block 


(17) 


3.2  Description  of  the  Monte  Carlo  Simulation 


The  objective  of  the  Monte  Carlo  simulation  was  to  grow  small 
initial  flaws  (0.005  inch  £ a^  <_  0.07  inch)  to  final  sizes  of  the 
order  of  the  critical  crack  size  (0.122  inch).  This  objective  was 
achieved  for  initial  flaws  in  the  range  0.03  <_  a^  £ 0.07  inch. 
However,  the  amount  of  growth  observed  for  very  small  flaws  (0.005 
£ ap  £ 0.01  inch)  was  insignificant  under  the  given  load  spectrum. 

Approximately  200  replicate  runs  were  made  for  each  a^  value. 
The  analog  output  of  crack  size  as  a function  of  time  was  recorded 
on  a strip  chart  recorder.  Subsequently,  templates  were  used  to 
read  the  crack  size  at  four  points  on  each  run,  corresponding  to: 

N = 8,000,  10,000,  13,300,  20,000  events  (18) 

The  crack  size  data  were  then  keypunched  and  processed  on  a digital 
computer  to  produce  observed  data  histograms  and  to  estimate  param- 
eters for  various  trial  distributions.  Data  histograms  for  final 
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crack  sizes  grown  from  0.03  <_  0.07  inch  are  shown  in  Figs.  15 

through  33. 

The  three-parameter  Weibull  distribution: 


, oC-l  _/X-X<A 

C,  V *(*-**)  „ 1 /*  1 

f60=  -J7T-  e 


fx-XoX 


oL 


Rx)=l-e  . * *>o  x^x0  (19) 


was  used  to  correlate  the  results.  The  initial  flaw  size  is  a 


logical  choice  for  the  lower-bound  parameter  xq  in  the  present 


case.  Hence,  the  model  becomes 


(20) 


where  a, 3 may  be  functions  of  a^  and  v,  and  where 

v = N/20,000  (21) 

defines  the  fraction  of  service.  The  crack  size  data  thus  corre- 
spond to  (see  Eq.  18): 

v = 0.4,  0.5,  0.667,  1 (22) 

Note  that  v = 1 is  equivalent  to  3,560  flight  hours. 

The  following  maximum  likelihood  estimator  formulas  for  a, 8 
may  be  derived  easily  for  the  three-parameter  Weibull  family  under 
the  assumption  that  the  lower-bound  parameter  is  known: 


A 

o(  = 


Nl  ( Xj : -X0) 

sj=1  4 


N X (X;-X0)*A t (*;-Xo) 


i4 


X (x.-xj 


(23) 


7 JU  (x:-X0) 

r'  J 


A 


(24) 
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where  x^ , X2,...,xN  are  the  observed  data.  Equation  23  may  be 
solved  for  a by  iteration,  after  which  3 is  obtained  explicitly 
from  Eq.  24.  Tables  3 and  4 summarize  the  estimates  thus  obtained 
for  a and  3,  respectively,  by  digital  processing  of  the  raw  data. 
The  case  of  a^  = 0.06  inch  (60  mils)  and  v = 1 is  an  exception,  the 

/V 

analog  data  having  been  accidentally  lost.  The  scale  parameter  3 
for  this  case  was  estimated  subsequently  from  the  empirical 
correlation  discussed  below. 


In  seeking  to  correlate  the  results  for  a and  3,  one  must  make 
some  approximate  hypothesis  which  can  be  followed  analytically  to  a 
conclusion  which  will  provide  a guide.  One  possible  hypothesis  is 
to  ignore  temporarily  the  random  nature  of  A a,  allowing  analytical 


integration  of  Eq.  10  to  the  form: 


for  Np  = 4 . If  this  hypothesis  is  combined  with  the  further 
assumption  that 


Pr{Final  Crack  Size  £ a^}^Pr{Load  £ Aa} 
and  Aa  is  modeled  with  a Weibull  distribution. 


p {&&)  = 1 - e 


there  results: 

| (a;.  |a0)>)  ^ F (txr)=  1-  exp* 


A 


(28) 
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Finally,  assuming  that  afa-j_'Vai  and  replacing  C^B  N by  Cv , where  C 
is  an  appropriate  constant,  leads  to: 


In  other  words,  the  shape  parameter  for  the  final  crack  size  model 
ought  to  be  about  one-fourth  of  the  estimated  loading  shape  param- 
eter, and  the  final  crack  size  scale  parameter  should  correlate 
2 

with  va^  . One  would  expect  this  chain  of  argument  to  hold  rigor- 
ously for  constant-amplitude  fatigue  with  randomly  chosen  amplitude, 
and  one  might  apply  it  approximately  if  most  of  the  crack  growth 
were  confined  to  the  last  few  load  cycles.  However,  its  utility 
in  the  present  case  is  limited  to  the  function  of  a general  guide. 

a 

Since  the  foregoing  argument  indicates  that  a constant  a is 
expected  for  final  crack  size  when  the  load  spectrum  is  kept  con- 

A 

stant,  mean  and  standard  deviation  were  estimated  for  a from  the 
data  in  Table  3,  with  the  results: 

= 1, 4?'  = O.ZZ*)  (30) 

These  results  may  be  compared  with  the  load  distribution  by  means 
of  the  well-known  relationship: 

* t/-/V  j * > i (31) 


Applying  Eq.  31  to  Eq.  9,  we  find: 


A « 14.25/8.25  = 1.73  (32) 

A 

Comparison  of  Eqs.  30  and  32  shows  that  the  estimated  mean  for  a 

A 

is  less  than  but  within  one  standard  deviation  of  A.  We  may 
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therefore  conclude  tentatively  that  the  variability  of  the  service 
load  distribution  is  impressed  more  or  less  directly  on  the  final 
crack  size  distribution. 

A 

The  final  crack  size  scale  parameters  B in  Table  4 were  sub- 


jected to  separate  regression  analyses  with  respect  to  v and  a^. 


These  analyses  indicated  that: 

ft  (vs-i') 


i.Z 


(33) 


The  final  correlation  for  B is  illustrated  in  Fig.  34 • The  correla- 
tion line  shown  in  the  figure  was  fit  by  eye  and  gives: 

js  - 2,l4*tO  / (34) 

with  a^  in  units  of  mils. 

The  foregoing  results  show  that  the  approximate  hypothesis 
based  on  a fully  integrated  growth  equation  does  indeed  provide 

/v 

only  a general  guide.  The  correlation  for  B is  in  general  agree- 
ment with  the  hypothetical  model.  However,  much  less  broadening 
of  the  final  crack  size  distributions,  relative  to  the  load  distri- 
bution, is  observed  in  comparison  with  the  prediction  of  the 
hypothetical  model.  For  the  purpose  of  further  analysis,  the  final 
crack  size  distribution  may  be  taken  as: 

^ r r % - et-  ~i 

* 1-  exP  j-  'l*.-s+*fo*(V9?)tA 

However,  no  special  significance  should  be  attached  to  Eq.  35.  It 
is  merely  the  result  of  an  empirical  correlation,  and  is  useful  only 
for  interpolation  within  and  near  the  ranges  of  initial  flaw  size  a^ 
and  service  fraction  V over  which  the  Monte  Carlo  experiments  were 


conducted. 
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3 . 3 Risk  Sensitivity  Assessment 


The  risk  for  a crack  to  exceed  critical  size  due  to  the  com- 
bined effects  of  random  loading  and  a random  initial  flaw  popula- 
tion may  be  assessed  by  combining  the  results  of  the  Monte  Carlo 
experiment  with  a description  of  the  initial  flaw  distribution. 
Since  it  is  logical  to  assume  that  the  random  variables  a^  and  A o 
are  uncorrelated,  the  following  simple  procedure  may  be  used  to 
compute  mathematical  risks: 

^ jp  ( ^cr  l V ^ “ fo  I ^ *i [’  (36) 

where  represents  the  risk  of  exceeding  acr,  given  the  service 
fraction  v,  and  where: 

} <37> 

f = Pr  f fUio  ^ ] (38> 

Approximate  risk  assessments  may  be  made  by  computing  Rf(acrlaj[/  v) 
at  a few  equally-spaced  a^  values  and  using  the  trapezoid  rule: 

Rifacr  !■*)  * C K (a,r|  j &*,*)  ■£.  (/  A*)  (39) 

0 

together  with  the  assumption  that  R^(a  |0,  v)  = 0. 

The  initial  flaw  distribution  referred  to  in  Eqs . 36,  38,  and 
39  must  be  defined  to  correspond  to  the  flaw  population  which  may 
escape  into  the  structure  at  the  beginning  of  service.  Models  for 
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this  population  may  be  constructed  by  combining  data  on  the  dis- 
tribution of  flaws  created  by  materials  processing  and  machining 
with  data  on  inspection  reliability.  Let  fm^g(a^)  be  the  prob- 
ability density  function  which  characterizes  the  "as-manufactured" 
flaw  population,  and  suppose  that  nondestructive  inspection  (NDI) 
data  are  available  in  the  form  of  a detection  probability  curve: 


Pd  (a)  = Pr{lnitial  flaws  >_  a are  detected  by  NDI} 


(40) 


(41) 


Some  typical  detection  probability  curves  are  shown  in  Fig.  35. 
Since  the  events  of  flaw  creation  and  detection  or  escape  may 
reasonably  be  assumed  to  be  uncorrelated,  the  initial  flaw  dis- 
tribution in  service  may  be  computed  from: 

p > = U-PjK)] 

The  effect  of  Ec.  41  is  illustrated  scehmatically  in  Fig.  36,  while 
Fig.  37  outlines  the  procedure  for  computation  of  the  mathematical 
risks  associated  with  exceedance  of  critical  crack  size.  Since 
there  are  commonly  no  data  available  for  estimation  of  (a) , 

good  statistical  practice  dictates  that  the  creation  of  flaws  of 
any  size  be  assumed  equally  likely  up  to  some  maximum  size  so  large 
that  the  risk  of  escape  from  detection  is  negligible.  In  this  case, 
Eq.  41  is  replaced  with: 

-p  (a)  = Va»i/x  ')  3-  £ 

"*3  (42) 


■f: 


(try  r^y 


l - P<)  M 

Ja 


(43) 
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An  example  risk  sensitivity  assessment  has  been  carried  out 
using  the  critical  crack  size  exceedance  function  R^(a^|a^,v) 
obtained  from  the  Monte  Carlo  simulation  in  combination  with  several 
possible  assumptions  about  the  initial  flaw  distribution  and  inspec- 
tion procedures.  Table  5 summarizes  the  critical  crack  exceedance 
risks  for  a matrix  of  service  fractions  v and  assumed  initial  crack 
sizes  a^.  The  results  in  this  table  include  extrapolation  outside 
the  a^  range  over  which  the  simulation  was  conducted.  In  Table  6 , 
the  details  of  the  computation  summarized  by  Eq.  39  are  shown,  using 
the  equal-likelihood  hypothesis  for  "as-manufactured"  initial  flaws 
(Eqs.  42  and  43)  and  the  "Mechanized  Probe"  detection  probability 
curve  from  Fig.  35.  The  summations  in  this  table  are  risks  associated 
with  equal-likelihood  as-manufactured  flaws,  mechanized-probe  NDI , 
and  a hypothetical  proof  test  which  positively  rejects  any  initial 
flaw  larger  than  100  mils.  To  obtain  the  corresponding  risks  for 
NDI  only,  one  need  simply  observe  that  Rf(acrlajyv)  ~ 1 for  aj_  > 100 
mils.  Hence,  the  case  of  NDI  only  may  be  computed  by  adding  the 
quantity 

1 - Fi(100  mils)  = 0.0428  (44) 

to  the  summations  in  Table  6 . 

It  is  interesting  to  compare  the  mathematical  risks  described 
above  with  a parallel  assessment  in  which  sensitivity  to  the  initial 
flaw  population  is  retained,  while  random  loading  effects  are  ignored. 
Since  the  random  process  adopted  for  the  above  examples  was  simply  a 
Gaussian-distributed  stress  range,  one  possible  measure  for  the 
corresponding  deterministic  loading  process  is  constant-amplitude 
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fatigue  at  the  spectrum  mean  stress  range  (14.25  ksi) . This  is 
referred  to  as  "central  loading"  for  convenience.  In  this  case, 
the  Paris  equation  may  be  integrated  to: 

-1 

(45) 

where  a^  and  a.  are  in  units  of  inches.  The  results  of  this  calcu- 
f l 

lation  are  summarized  in  Table  7.  Critical  crack  exceedance  risks 

for  the  case  of  central  loading  may  be  computed  as  follows.  For  a 

given  service  fraction  v,  determine  by  interpolation  from  Table  7 

the  a.  value  for  which  a,  = a = 122  mils.  The  risk  is  then 
1cr  f cr 

given  by  1 - F . (a . ) . 

cr 

Finally,  three  additional  cases  have  been  included,  in  which 
some  hypothetical  assumptions  have  been  made  about  the  as-manufactured 
flaw  population  and  for  which  UDI  without  a proof  test  has  been 
assumed.  Hypothetical  Gaussian  distributions: 

= } '46> 
were  assumed,  with  a in  units  of  mils,  and  with  a = 10,  15,  20  mils. 
These  distributions  thus  have  coefficients  of  variation  of  0.2,  0.3, 
0.4  respectively.  Table  8 summarizes  the  as-manufactured  probability 
densities  and  the  initial  flaw  probability  densities  computed  from 
Eq.  41.  Combined  risk  calculations  have  been  carried  out  for  each  of 
these  distributions,  as  shown  in  Tables  9,  10,  and  11. 

Table  12  gives  a comparative  summary  of  the  mathematical  risks 
computed  for  each  of  the  cases  described  above.  The  first  two  lines 
in  this  table  assume  an  equal-likelihood  manufactured  population  and 
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account  for  random  loading.  The  third  and  fourth  lines  correspond 
to  deterministic-load  risk  assessment  with  the  same  initial  popula- 
tions. The  last  three  lines  account  for  randomness  in  the  loading 
and  the  initial  flaw  population,  and  are  intended  to  illustrate  the 
sensitivity  of  risk  to  dispersion  in  the  initial  population. 

The  mathematical  probabilities  in  Table  12  present  quite  different 
pictures  of  the  risk  of  achieving  a critical  crack  size.  At  one 
extreme,  two  models  predict  small  or  no  risk  at  short  service  times, 
with  the  risk  function  rising  rapidly  as  the  service  approaches  one 
life.  These  two  extreme  models  are  the  central  load  calculation 
with  equal-likelihood  initial  population  filtered  by  NDI  and  proof 
test,  and  the  random  load  calculation  with  Gaussian  initial  popula- 
tion (CV  = 0.2)  filtered  by  NDI.  At  the  other  extreme,  two  models 
predict  rather  high  risk  over  the  entire  service  period  considered, 
with  the  risk  function  relatively  insensitive  to  service  time.  These 
are  the  random  load  calculation  with  equal-likelihood  initial  popula- 
tion filtered  by  NDI,  and  the  corresponding  central  load  calculation. 
The  other  cases  are  seen  to  give  risk  functions  which  lie  between 
these  extremes. 

It  is  apparent  from  these  results  that  the  procedures  which 
might  be  adopted  to  minimize  risks  for  service  times  less  than  one 
life  depend  strongly  upon  which  model  is  accepted  as  a true  descrip- 
tion of  nature.  Thus,  if  the  equal-likelihood  population  is  assumed 
for  manufactured  flaws  and  random  load  is  accounted  for,  one  is  led 
to  strategies  of  either  a proof  test  to  reject  initial  flaws  much 
smaller  than  100  mils  or  development  of  some  NDI  method  which  will 
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give  better  detection  probabilities.  On  the  other  hand,  the  given 
NDI  method  supplemented  by  the  100-mil  proof  test  appears  to  be  a 
.good  strategy  if  the  effects  of  random  loading  are  ignored.  Of 
course,  the  equal-likelihood  manufactured  flaw  population  is  a very 
conservative  assumption,  even  though  accepted  statistical  practice 
demands  it  in  the  absence  of  experimental  data.  The  examples  with 
Gaussian  manufactured  populations  demonstrate  that  it  may  be  worth- 
while to  obtain  enough  data  to  determine  at  least  approximately  the 
coefficient  of  variation  for  as-manufactured  flaw  populations. 
Finally,  it  is  interesting  to  compare  these  latter  cases  directly 
with  the  data  in  Table  5 . Each  line  in  Table  5 may  be  treated  as 
a risk  function  which  corresponds  to  an  assumption  that  all  initial 
flaws  have  a common  deterministic  size,  a view  which  is  convenient 
for  design  calculations.  The  question  is  whether  it  is  possible  to 
choose  a flaw  size  which  leads  to  a risk  function  with  behavior 
similar  to  the  risk  functions  based  on  the  entire  flaw  population. 

In  the  present  case,  the  answer  is  affirmative  if  a^  * 70  mils  is 
chosen  for  short  times  and  a^  - 55  mils  is  chosen  for  service  periods 
approaching  one  life. 
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Section  4 


DISCUSSION  AND  CONCLUSIONS 


Results  have  been  presented  for  a fatigue  crack  growth  risk 
analysis  for  which  a Gaussian  stress  range  spectrum  was  assumed, 
and  for  which  a critical  crack  size  for  fracture  was  estimated 
approximately  as  0.122  inch  in  7075-T6  aluminum  alloy.  The  stress 
range  spectrum  was  found  to  be  capable  of  growing  cracks  from 
initial  sizes  between  0.03  and  0.07  inch  to  the  critical  crack 
size  after  20,000  stress  range  events,  a total  which  corresponds 
to  approximately  3,600  flight  hours. 

The  objective  of  the  analysis  was  to  generate  a data  base  for 
final  crack  size  at  several  fractions  of  service,  from  which  the 
risks  of  structural  failure  could  be  assessed.  The  risk  study 
focussed  primarily  on  the  contribution  of  dispersion  in  the  load 
spectrum  to  the  dispersion  of  final  crack  sizes.  Average  service 
experience  for  a fleet,  when  translated  into  a final  crack  size, 
provides  only  limited  information  about  fleet  safety.  The  calcula- 
tion of  failure  rates  requires  crack  growth  computations  for  non- 
representative as  well  as  average  service  history  within  a given 
mission  spectrum.  Since  no  analytical  derivation  can  be  given  for 
the  distribution  of  final  crack  sizes  as  a function  of  the  stress 
range  distribution,  a Monte  Carlo  simulation  is  required,  in  which 
a crack  growth  rate  equation  is  integrated  several  hundred  times 
from  the  same  initial  crack  size.  A technique  of  analog  computer 
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simulation  was  developed  and  verified  for  this  purpose.  The  analog 
method  was  shown  to  be  cost-effective  for  extended  Monte  Carlo 
simulations,  in  comparison  to  digital  computation.  Also,  the 
stress  range  spectrum  was  shown  to  be  easily  representable  by  input 
of  a filtered  signal  from  a random  noise  generator.  The  random  noise 
signal  automatically  provides  a series  of  load  histories  with  con- 
trolled parameters  which  fit  the  given  mission  spectrum  (stress 
range  mean  and  standard  deviation)  as  the  analog  repeatedly  integrates 
the  crack  growth  rate  equation.  With  the  analog  set  in  repetitive 
operation,  approximately  30  minutes  were  required  to  generate  200 
curves  of  crack  size  versus  equivalent  flight  hours  for  each  of  11 
assumed  initial  crack  sizes.  Four  data  points  were  sampled  from 
each  curve  (0.4,  0.5,  2/3  and  full  service  period)  and  were  processed 
digitally  to  provide  histograms  for  final  crack  size,  with  time  and 
initial  size  as  parameters.  The  raw  data  histograms  for  five  initial 
sizes  are  presented  in  this  report. 

The  three-parameter  Weibull  distribution  was  chosen  as  a logical 

model  for  the  probability  distribution  for  final  crack  size.  With 

the  lower-bound  location  parameter  assumed  equal  to  the  initial  crack 

size  for  each  histogram,  the  analog  data  were  used  to  compute  maximum 

likelihood  estimates  for  the  scale  and  shape  parameters.  A regression 

analysis  of  the  scale  parameters  resulted  in  good  empirical  correla- 
2 12 

tions  with  (va^)  ’ , where  v represents  service  fraction  and  a^  is 

the  initial  crack  size.  These  results  agreed  generally  with  an 

2 

approximate  estimate  which  predicted  va^  for  the  correlation.  The 
approximate  estimate  also  predicted  a final  crack  size  dispersion 
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of  four  times  the  load  spectrum  dispersion.  However,  the  final 
crack  size  shape  parameters  estimated  from  the  analog  data  indicated 
that  the  dispersion  of  final  crack  size  is  only  slightly  greater 
than  the  load  dispersion.  Shape  parameters  for  final  crack  size 
were  therefore  defined  accordingly,  and  were  used  together  with 
the  empirical  correlations  for  scale  parameter  to  extrapolate  the 
data  base  to  a few  near-neighboring  initial  crack  sizes  which  had 
not  been  included  in  the  simulation. 

The  final  crack  size  data  base  was  used  to  compute  the  mathe- 
matical risks  of  exceeding  the  critical  crack  size,  as  a function 
of  service  time,  for  several  sets  of  assumptions  about  initial 
flaw  populations  and  manufacturing  quality-control  strategies. 

For  the  parameters  considered  in  the  present  study,  the  results 
indicated  that  risk  functions  based  on  deterministic  assumed  initial 
flaw  sizes  between  55  and  70  mils  were  similar  to  the  risk  functions 
based  on  assumed  flaw  populations  with  a mean  size  of  50  mils  and 
coefficients  of  variation  between  0.2  and  0.4  which  were  subjected 
to  filtering  by  nondestructive  inspection.  Risk  functions  were  also 
computed  for  cases  which  considered  a flaw  population  but  neglected 
random  load  effects.  These  risk  functions  were  observed  to  behave 
quite  differently  from  those  calculated  for  corresponding  cases  which 
included  the  effects  of  random  loading.  This  comparison  may  have 
a significant  bearing  on  risk  assessments  of  aircraft  which  are 
fatigued  primarily  by  gust  loading,  but  must  be  considered  much 
less  important  for  aircraft  which  are  fatigued  primarily  by  quasi- 
deterministic  maneuver  loads. 
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There  appears  to  be  no  reason  why  the  analog  method  cannot  be 
extended  to  more  sophisticated  simulations,  in  which  the  effects 
of  stress  ratio,  retardation,  multiple  sources  of  random  loading, 
ground-air-ground  cycles,  etc.,  are  included.  The  feasibility  of 
such  extensions  depends,  of  course,  upon  the  availability  of  a 
sufficient  amount  of  analog  hardware,  and  upon  calibration  by 
parallel  digital  analysis  to  verify  computational  accuracy.  Given 
these  conditions,  the  most  useful  application  of  the  analog  approach 
to  crack  growth  simulation  appears  to  be  rapid  assessment  of  the 
sensitivity  of  fleet  risk  functions  to  changes  in  load  spectra. 
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TABLE  1 


COMMON  LINEAR  ANALOG  COMPONENTS 


SIMPLE  CIRCUITS  USING  AMPLIFIERS  AND  POTENTIOMETERS 


Courtesy  of  EAI-TR-48  Reference  Guide 


*2e 


Courtesy  of  EAI-TR-48  Reference  Guide 
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TABLE  2 


COMMON  NONLINEAR  ANALOG  COMPONENTS 


QUARTER- SQUARE  MULTIPLIER  CIRCUITS 


PATCHING  DIAGRAM 


Courtesy  of  EAI-TR-48  Reference  Guide 
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Courtesy  of  EAI-TR-48 


CIRCUITS 


Guide 


Courtesy  of  EAI-TR-48  Reference  Guide 
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TABLE  3 


a . 

l 

(mils) 


WEIBULL  SHAPE  PARAMETER  a FOR 
FINAL  CRACK  SIZE  DISTRIBUTION 


N = 8 Kc 
v = 0.4 


10  Kc  13.3  Kc 

0.5  0.667 


20  Kc 
1 


30 

1.61 

1.76 

1.99 

2.16 

40 

1.36 

1.49 

1.63 

1.71 

50 

1.28 

1.33 

1.49 

1.70 

60 

1.12 

1.20 

1.33 

* 

70 

1.41 

1.45 

1.15 

1.00 

*Data 

for  a.  = 
i 

60  mils,  N = 

20  Kc  were 

lost. 
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TABLE  4 


/\ 

WE I BULL  SCALE  PARAMETER  3 (mils)  FOR 
FINAL  CRACK  SIZE  DISTRIBUTION 


a . 
i 

N = 8 Kc 

10  Kc 

13.3  Kc 

20  Kc 

(mils) 

v = 0.4 

in 

• 

o 

0.667 

1 

30 

2.71 

3.38 

4.58 

7.55 

40 

4.15 

5.41 

8.01 

13.42 

50 

10.17 

12.72 

19.12 

28.63 

60 

12.42 

15.84 

23.16 

43.5* 

70 

18.72 

23.51 

37.43 

71.35 

♦Estimated  from  empirical  correlation 
g=2.34xl0“3(va?)1*2 
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TABLE  5 


RISK  R, (a  a . , v ) FOR  a = 122  MILS 

31  Ci  X C J- 


Service  Fraction  v 

(mils) 


0.4 

0.5 

2/3 

1 

20 

-391 

1.00x10 

l.OOxlO-261 

3 . 9 8xl0~157 

1.00x10 

30 

_ . , ,,,-80 
3 . 16x10 

1.00x10  33 

7.08xl0_33 

3 . 80x10 

40 

-24 

1.26x10 

7 . 94xl0~17 

1.82xl0-1° 

2 . 02x10 

50 

1. 26xl0_9 

8 . 21xl0~  7 

2.63xl0-4 

1.74x10 

60 

1. 74xl0_4 

2.88xl0~3 

2.89xl0~2 

1.84x10 

70 

2.14xl0-2 

7.51xl0_2 

2.14xl0-1 

4.65x10 

80 

1 . 74xl0_1 

3. 07xl0_1 

4. 91xl0_1 

7.08x10 

90 

4.64x10  3 

5 . 95xl0_1 

7. 33xl0-1 

8.60x10 

100 

7. 31xl0_1 

8 . 16xl0_1 

8. 85xl0-1 

9.41x10 
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TABLE  6 


COMBINED  RISK  FOR  a±  TO  acr  = 122  MILS 


Products  Rf (acr|ai,v)fi(ai) Aa  at  Fraction  v 

ai  fi(ai)Aa* 


(mils) 

0.4 

0.5 

2/3 

1 

20 

.1602 

1.60xl0-392 

1.60x10  262  6.38x10  158 

1.60X10-77 

30 

.1010 

3 . 19xl0~81 

-54 

1.01x10 

7 . 15xl0~34 

3 . 84xl0~17 

40 

.0593 

7.47xl0“26 

4 . 7 0xl0~18 

1.08xl0-11 

1.20xl0-6 

50 

.0321 

4.04xl0_11 

2.64xl0~8 

8 . 44xl0~6 

5.59xlO"4 

60 

.0192 

3.34xl0~6 

5.53xl0-5 

5. 55xl0~4 

3. 53xl0-3 

70 

.0128 

2.74x10  4 

9 . 61xl0~4 

2 . 74xl0~3 

5.95xl0"3 

80 

.0112 

1 . 95xl0~3 

3.44xl0-3 

5 . 50xl0~3 

7 . 93xl0~3 

90 

.00961 

4.45x10  3 

5 . 71xl0~3 

7. 04xl0~3 

8.26xl0-3 

100 

.00801 

5.85x10  3 

6 . 54xl0-3 

7.09x10  3 

7 . 54xl0-3 

Z 

1.25xl0“2 

1 . 67xl0~2 

2 . 29xl0~2 

3 . 38xl0~2 

Computed 

from  " 

Mechanized  Probe 

" curve 

(Fig.  35)  . 
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TABLE  7 


CRACK  GROWTH  BEHAVIOR  UNDER  CENTRAL  LOADING 


a . 
1 

(mils) 

at  Service  Fraction  v 

(mils) 

0.4 

0.5 

2/3 

1 

20 

20.5 

20.7 

20.9 

21.4 

30 

31.2 

31.6 

32.0 

33.2 

40 

42.2 

42.8 

43.7 

46.0 

50 

53.4 

54.5 

56.0 

59.6 

60 

65.0 

65.5 

68.8 

74.4 

70 

76.9 

78.9 

82.3 

90.4 

80 

89.1 

92.0 

96.5 

108.0 

90 

102.0 

105.4 

111.6 

127.0 

100 

115.0 

119.2 

127.0 

110 

128.0 

133.5 
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TABLE  8 


HYPOTHETICAL  DISTRIBUTIONS  FOR 
INITIAL  FLAWS 


a . 

^mf  g 

(a . ) Aa 

l 

f . (a . ) Aa 

l l 

(mils) 

CM 

O 

II 

> 

a 

CV=0. 3 

o 

II 

> 

CJ 

CM 

O 

II 

> 

CJ 

CV=0. 3 

o 

II 

> 

CJ 

10 

— 

.0114 

. 0540 

— 

.0314 

.1063 

20 

.0044 

.0540 

.1295 

.0140 

.0985 

.1685 

30 

.0540 

.1828 

.2420 

.1110 

.2140 

.2020 

40 

. 2420 

.3256 

.3521 

.3200 

.2475 

.1912 

50 

.3989 

. 3989 

. 3989 

. 3600 

. 2075 

.1440 

60 

. 2420 

. 3256 

.3521 

.1640 

.1265 

.0977 

70 

.0540 

.1828 

.2420 

. 0291 

.0568 

.0538 

80 

.0044 

.0540 

.1295 

. 0020 

.0140 

.0239 

90 

1. 36xl0-5 

.0114 

.0540 

5 . 04x10 

6 .0024 

.00824 

100 

1.49xl0"7 

3 . 28xl0~6 

.0175 

4 . 50x10 

8 5.74xl0_7 

.00217 

110 

— 

— 

.0044 

— 

— 

.00044 

120 

— 

4 

.38xl0~5 

— 



3.60x10' 
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TABLE  9 


COMBINED  RISK  FOR  a.  TO  a (CV=0.2) 

X 


Products  Rf 

(a  la. fv) f . 

(a.)Aa  at  Fraction  v 

cr  1 i i 

a . 

f . (a . ) Aa 

1 

l l 

(mils) 

(Table  12) 

0.4 

0.5 

2/3 

1 

40 

. 3200 

4 . 03xl0~25 

2 . 54xl0~17 

5. 82xl0_11 

6.46x10  6 

50 

.3600 

4.54xl0_1° 

2 . 95xl0~7 

9.47xl0-5 

6.26xl0-3 

60 

.1640 

2. 86xl0~5 

4.72xl0_4 

4.74xl0-3 

3 . 02xl0~2 

70 

.0291 

6.24x10  4 

2 . 18xl0~3 

6.24xl0-3 

1.35xl0-2 

80 

. 0020 

3.48xl0-4 

6. 14xl0'4 

9.82xl0_4 

1.42xl0-3 

90 

5 . 04xl0~  6 

2. 34xl0-6 

3. OOxlO-6 

3 . 69xl0~6 

4.34xl0-6 

100 

4.50xl0_8 

3. 29xl0~8 

3. 68xl0-8 

3.98xl0-8 

4.24xl0-8 

E 

1. OOxlO-3 

3.27xl0-3 

1.21xl0-2 

5.14xl0-2 
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TABLE  10 


COMBINED  RISK  FOR  a . TO  a (CV=0.3) 

x cr 


Products  R,-  (a  a . , v ) f . (a  . ) Aa  at  Fraction  v 
f cr  x xx 


a . 

X 

f . (a . ) Aa 

X X 

(mils) 

(Table  12) 

40 

. 2475 

3 . 12xl0~25 

1.96xl0“17 

4.50xl0_11 

5.00xl0-6 

50 

.2075 

2. 62xl0-10 

1. 71xl0~7 

5.46xl0”5 

3.61xl0-3 

60 

.1265 

2. 20xl0_5 

3. 64xl0~4 

3. 66xl0-3 

2.33xl0-2 

70 

.0568 

1. 21xl0~3 

4 . 26xl0~3 

1.21xl0“2 

2. 64xl0-2 

80 

.0140 

2 . 44xl0~3 

4.30xl0“3 

6 . 86xl0~3 

9.90xl0-3 

90 

.0024 

l.llxlO-3 

1.43xl0-3 

1. 76xl0-3 

2.06xl0-3 

100 

5 . 74xl0~7 

4.20xl0“7 

4. 63xl0~7 

5 . 08xl0~7 

5 . 41xl0~7 

£ 

4.78xl0~3 

1 . 04xl0~2 

2 . 44xl0~2 

6 . 53xl0~2 

t 
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TABLE  11 


COMBINED  RISK  FOR  a . TO  a (CV=0.4) 

l cr 


Products  Rj.  (a  I a . , v ) f . (a  . ) Aa  at  Fraction  v 

f cr 1 1 i i 


a . 

l 

(mils) 

f . (a . ) Aa 

l l 

(Table  12) 

0.4 

0.5 

2/3 

1 

50 

.1440 

1.81xl0"10 

1 . 18xl0~7 

3.79xl0“5 

2 . 51xl0~3 

60 

.0977 

1.70xl0“5 

2. 81xl0~4 

2. 82xl0_3 

1. 80xl0~2 

70 

.0538 

1 . 15xl0_3 

4 . 04xl0~3 

1.15xl0_2 

2 . 50xl0~2 

80 

.0239 

4.16xl0"3 

7. 35xl0~3 

1 . 17xl0~2 

1 . 69xl0~2 

90 

.00824 

3. 82xl0~3 

4.90xl0~3 

6 . 04xl0_3 

7 . 09xl0~3 

100 

. 00217 

1. 59xl0~3 

1 . 77xl0~  3 

1.92xl0-3 

2 . 04xl0~3 

110 

.00044 

4 . 40xl0~4 

4. 40xl0~4 

4 . 4 0xl0_4 

4. 40xl0~4 

E 

1 . 12xl0~2 

1. 88xl0~2 

3. 45xl0-2 

7 . 20xl0~2 
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0 o.l  0.2  0.3 

CRACK  SIZE  a INCH 


Figure  1:  Weibull  Distribution  for  Crack  Size  as  a Function  of 

Constant-Amplitude  Fatigue 


Oo/ 

//o 


Figure  2:  Paris  Equation  Analog  Program  Diagram 
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~7~£%T  Ru\s  FdF  CoK  sfaajT  Check 


Figure  3:  Typical  Analog  Crack  Growth  Solution  for  Constant- 

Amplitude  Fatigue 
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J x/0 

((roe  cYcees 


Track-  Hold  ! cm  ~ 0 , Z m s /Doo  Hz 


Half  (Kick  - Half  Houd 

77*ack- Hold  l Cm  - 10  ms  loo  Hi 


Max  Hold 

Figure  5:  Samples  of  Track-Hold  Noise 
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Figure  7:  Analog  Program  Diagram  for  Paris  Equation  with 

Random  Loading 
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REPOP  Compute  77/ mb  ~ 2-0 0 nsec  / PfioBU?" 


Figure  10:  REPOP  Switch  Error 
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Figure  11:  Analog  Program  Diagram  for  Paris  Equation  with 
Random  Loading  and  Rescaling 
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O-iO  0.20 

CRACK  SIZE  a inch 


Figure  12:  Histogram  of  Crack  Size  after  10,000  Fatigue 
Cycles 
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0.20 


0 . 30 


CRACK  SIZE  a inch 


Figure  13:  Histogram  of  Crack  Size  after  15,000  Fatigue 
Cycles 
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NUMBER 
OF  DATA 


0.03  0.04 

CRACK  SIZE  a inch 

a = 0.03"  N = 8,000  FATIGUE  CYCLES 


Figure  15:  Crack  Size  Histogram  (a. 


.03,  N 


8,000) 
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CRACK  SIZE  a inch 

0.03"  N = 10,000  FATIGUE  CYCLES 


Figure  16:  Crack  Size  Histogram  (a.. 


.03,  N = 10,000) 
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0.03  0.04  0.05 

CRACK  SIZE  a inch 

0.03"  N - 13,300  FATIGUE  CYCLES 


Figure  17:  Crack  Size  Histogram  (a±  = .03,  N = 13,300) 
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NUMBER 
OF  DATA 


CRACK  SIZE  a inch 

aQ  = 0.03"  N = 20,000  FATIGUE  CYCLES 


Figure  18:  Crack  Size  Histogram  (a.  = .03,  N = 20,000^ 
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CRACK  SIZE  a inch 

aD  = 0.04"  N = 8,000  FATIGUE  CYCLES 


Figure  19:  Crack  Size  Histogram  (a ^ = .04,  N = 8,000) 
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a = 0.04"  N = 20,000  FATIGUE  CYCLES 
o 

Figure  22:  Crack  Size  Histogram  (ai  = .04,  N -•  20,000) 
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Figure  23:  Crack  Size  Histogram  (a^  = .05,  N = 8,000) 
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0.05  0.10  0.15 

CRACK  SIZE  a inch 

a = 0.05"  N = 13,300  FATIGUE  CYCLES 

o 


Figure  25:  Crack  Size  Histogram  (a. 


.05,  N = 13,300) 


CRACK  SIZE  a inch 

aQ  = 0.05"  N = 20,000  FATIGUE  CYCLES 


Figure  26:  Crack  Size  Histogram  (a^  - .05,  N - 20,000) 
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N = 8,000  FATIGUE  CYCLES 


Figure  27:  Crack  Size  Histogram  (a. 


.06,  N = 8,000) 


CRACK  SIZE  a inch 

aQ  = 0.06"  N = 13,300  FATIGUE  CYCLES 

Figure  29:  Crack  Size  Histogram  (a.  = .06,  N = 13,300) 
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aQ  = 0.07"  N = 8,000  FATIGUE  CYCLES 


Figure  30:  Crack  Size  Histogram  (a. 


.07,  N = 8,000) 


Figure  31:  Crack  Size  Histogram  (a.  = .07,  N - 10,000) 
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o 

Figure  32:  Crack.  Size  Histogram  (a.  = .07,  N = 13,300) 
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Figure  34:  Correlation  of  Final  Crack  Size  Distribution 
Parame ters 
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Figure  35:  Typical  Detection  Probability  Curves 


Final  Crack  Size,  a. 


Initial  Crack  Size,  a. 


Figure  37:  Schematic  for  Computation  of  Combined  Failure  Risk 
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